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Abstract 

Graphene plasmons provide a suitable alternative to noble-metal plasmons because 
they exhibit much larger confinement and relatively long propagation distances, with 
the advantage of being highly tunable via electrostatic gating. We report strong light- 
matter interaction assisted by graphene plasmons, and in particular, we predict un- 
precedented high decay rates of quantum emitters in the proximity of a carbon sheet, 
large vacuum Rabi splitting and Purcell factors, and extinction cross sections exceed- 
ing the geometrical area in graphene ribbons and nanometer-sized disks. Our results 
provide the basis for the emerging and potentially far-reaching field of graphene plas- 
monics, offering an ideal platform for cavity quantum electrodynamics and supporting 
the possibility of single-molecule, single-plasmon devices. 
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Surfaces plasmons (SPs), the electromagnetic waves coupled to charge excitations at 
the surface of a metal, are the pillar stones of applications as varied as ultrasensitive opti- 
cal biosensing/"^ photonic metamaterials,^ light harvesting, ^'^ optical nano-antennas,^ and 
quantum information processing. However, even noble metals, which are widely regarded 
as the best available plasmonic materials, are hardly tunable and exhibit large ohmic losses 
that limit their applicability to optical processing devices. 

In this context, doped graphene emerges as an alternative, unique two-dimensional plas- 
monic material that displays a wide range of extraordinary properties. This atomically 
thick sheet of carbon is generating tremendous interest due to its superior electronic and 
mechanical properties, which originate in part from its charge carriers of zero effective 
mass (the so-called Dirac fermions"^^) that can travel for micrometers without scattering, even 
at room temperature.^^ Furthermore, rapid progress in growth and transfer techniques have 
sparked expectations for large-scale production of graphene-based devices and a wide range 
of potential applications such as high-frequency nanoelectronics, nanomechanics, transparent 
electrodes, and composite materials. 

Recently, graphene has also been recognized as a versatile optical material for novel 
photonic and optoelectronic applications,^^ such as solar cells, photodetectors,^^ light 
emitting devices, ultrafast lasers, optical sensing, and metamaterials.^^ The outstanding 
potential of this atomic monolayer is emphasized by its remarkably high absorption ^^'^^ 
^ ira ^ 2.3%, where a — /he ^ 1/137 is the fine-structure constant. Moreover, the linear 
dispersion of the Dirac fermions enables broadband applications, in which electric gating can 
be used to induce dramatic changes in the optical properties. 

All of these photonic and optoelectronic applications rely on the interaction of propa- 
gating far-field photons with graphene. Additionally, SPs bound to the surface of doped 
graphene exhibit a number of favorable properties that make graphene an attractive alter- 
native to traditional metal plasmonics. In particular, graphene plasmons are confined to 
volumes ^10^ times smaller than the diffraction limit, thus facilitating strong light-matter 
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interactions. Furthermore, dramatic tuning of the plasmon spectrum is possible through 
electrical or chemical modification of the charge carrier densities, ^^'^^ and a Fermi energy 
as high as =1-2 eV has been recently realized. ^^'^-^ And lastly, the electronic structure 
of graphene and the ability to fabricate large, highly crystalline samples should give rise 
to SP lifetimes reaching hundreds of optical cycles, thereby circumventing one of the major 
bottlenecks facing noble-metal plasmonics. 

Here, we show that these properties can be used to tailor extremely strong light-matter 
interactions at the quantum level. In particular, we consider the interaction between a single 
quantum emitter and single SPs in graphene, and show that the extreme mode confine- 
ment yields ultrafast and efficient decay of the emitter into single SPs of a proximal, doped 
graphene sheet (see Fig. la). More precisely, we analyze confinement in 2D (homogeneous 
graphene), ID (nanoribbon), and OD (nanodisk) geometries. We find an increased degree of 
field enhancement and interaction strengths with reduced dimensionality, ultimately yielding 
in OD structures decay rates exceeding the natural decay rate by six orders of magnitude. 
Consequently, graphene opens up a novel route to quantum plasmonics and quantum devices 
that have so far been difficult to achieve with conventional metal plasmonics. Beyond the 
controlled enhancement and channeling of emission, graphene should also assist the explo- 
ration of fundamentally new regimes of quantum plasmonic interactions at the nanoscale. 
Specifically, we predict observable vacuum Rabi splittings in our proposed nanostructures, 
enabling a single SP to be emitted and then re-absorbed. Finally, we show that these 
structures have resonant extinction cross-sections greatly exceeding their geometrical cross- 
sections, despite the small volume occupied by this thin material, thus rendering the effects 
observable in practice and paving the way to advanced optoelectronic applications in which 
photon absorption is dramatically enhanced. 
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Optical response of graphene 

The photonic properties of this material can be fully traced back to its in-plane conduc- 
tivity (j(fc||,CL;), which is in general a function of parallel wave vector fcy and frequency uo. 
This quantity is mainly controlled by electron-hole pair (e-h) excitations that can be di- 
vided into intraband and interband contributions (see Fig. lb). Within the random-phase 
approximation ^^'^^ (RPA), the conductivity of graphene in the local hmit {k\\ 0) reduces 
to 



a{uj) = — : — ^ + — - 



9{huj-2EF) + -log 

TT 



(1) 



huj + 2EF 

(we take Ep > 0, see Fig. Ic). The first term of Eq. (Eq. (1)) describes a Drude model 
response for intraband processes, conveniently corrected to include a finite relaxation time 
r, for which we use a conservative value (r ^ 10""^^ s at Ejp — 0.1 eV~"^) extracted from the 
measured, impurity-limited DC mobility ■'^^'■'^^ (see Appendix). The remaining terms arise 
from interband transitions, which produce significant losses at energies near and above 2Ef 
(i.e., Rejcr} = resulting in the well-established 2.3% absorption). Equation (Eq. (1)) 

is valid at zero temperature, but we actually employ a finite-temperature extension^^ at 
T = 300 K (see Appendix), in which the step function is smeared out by thermal effects (Fig. 
Ic). The local hmit (Eq. (Eq. (1))) produces decay rates in reasonable agreement with the 
nonlocal RPA (see Appendix), except in the region near the interband transition onset or 
for distances below vp/uj. 



Plasmons in homogeneous graphene: Extraordinary confinement 

For sufficiently high doping [Ep > uj) graphene can sustain p-polarized SPs propagating along 
the sheet with wave vector /cgp ~ z(e+ l)cj/47rcr and electric field profile E ^ exp[/csp(zx — \z\)]. 
These electrostatic expressions, valid for hu ^ aEp (see Appendix), reduce upon insertion 
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of the Drude formula (first term of Eq. (Eq. (1))) to 

fcsp ^ {h^/Ae^EF) (e + 1) + i/t) , (2) 

clearly showing a quadratic dependence of fcgp on cj, which is characteristic of 2D electron 
gases (see Fig. le). 

The remarkable degree of confinement provided by the graphene is clear from the ratio of 
SP to free-space-light wavelengths Agp/Ao ^ [Aa / {e+l)\{EF / huo) derived from Eq. (Eq. (2)). 
In addition, the out-of-plane wave vector ^ zfcgp indicates an equally tight confinement to 
dimensions ^ Asp/27r in the transverse direction z. 

Interestingly, the in-plane propagation distance (1/e decay in amplitude), given by l/Imj/cgp}, 
reaches values well above 100 SP wavelengths (see inset of Fig. le) and drops rapidly at high 
energies when the plasmon has sufiicient energy to generate e-h's and the dispersion relation 
enters the intraband region (see Fig. lb). Figure le shows clearly that for increasing Bp the 
plasmons become narrower because the damping rate 1/r decreases relative to the photon 
frequency uj (Fig. le). 

Strong SP-emitter coupling in graphene 

Our first realization is that when an emitter such as an excited molecule or a quantum dot is 
placed close to doped homogeneous graphene, the emission rate is enhanced by 1-5 orders of 
magnitude and the emitted energy is mainly converted into a plasmon in the carbon sheet, 
as illustrated in the near-electric-field-amplitude plot of Fig. la. 

The decay rate F is proportional to the strength of the coupling between the transition 
dipole matrix element d and the electromagnetic modes acting on it, including the plasmon. 
This can be related to the electric field induced by the dipole on itself E^^^ (i.e., the field 
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reflected by the graphene) as^^ 

r = ro + |lm{d*-Ei'^d|^ 

I L 

where Fq = 4A;Q|dp/3/i is the free-space decay rate. For homogeneous graphene, the induced 
field is related to the Fresnel coefficients. This yields an exact relation for F as an integral 
over parallel wave vector contributions^^ (see Appendix), which we use in the calculations 
presented in Figs. 1,2. However, it is instructive to explore the electrostatic limit, which is 
accurate for small distances compared to the emission wavelength: 



r ^ Fo + I f |d|| |2 + 2|d^|2) r kl dku Im I j^—— 



e-2^ll^, (4) 



\\a/uj 



where z is the emitter- graphene separation and || (_L) denotes components parallel (perpen- 
dicular) to the graphene. The exponential in the above integral effectively suppresses the 
contribution of wave vectors /cy ^ 1/z. 

The spectral dependence of the decay rate is represented in Fig. Id (sohd curves) for 
various values of when the emitter is placed 10 nm away from graphene supported on 
silica (e = 2). The rate is peaked at a photon energy below Ep^ before dropping dramatically 
and finally converging to a common £^i?-independent value at energies above 2Ef. This 
behavior can be understood as follows. When the SP mode is well defined, the integral in 
Eq. (Eq. (4)) separates into two distinct contributions, a sharp pole associated with emission 
into the SPs and a broad background out to wave vectors fcy ^ 1/z associated with emission 
into lossy channels. The pole contribution yields the SP emission rate 



F, 



This SP contribution (dashed curves in Fig. Id) dominates the decay below Ep and is 
responsible for the emission maximum in that region. For distances z <C Agp, the decay rate 
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is enhanced by a factor [37r//2(e + 1)] (Ao/Agp)^ relative to the free-space rate, where Aq 
Asp) is the hght wavelength and / = 1 (/ = 2) for parallel (perpendicular) polarization. 
This explains the high values obtained for T/Tq in Fig. Id, which can be attributed to the 
smaller mode volume ((Ao/Agp)^ factor) and the reduced group velocity (additional Ao/Agp 
factor) . 

The noted spectral dip in the decay corresponds to the onset of interband transitions at 
hcu — 2Ef and decay into lossy channels as the dominant emission mechanism (see Fig. lb,c). 
At energies above this dip, the response is dominated by interband transitions, whereas intra- 
band and SP excitations become unimportant, and the rate follows a common profile similar 
to undoped graphene [Ep = 0). In fact, undoped graphene exhibits a novel phenomenon of 
strong quenching induced by the high conductivity of the carbon sheet. 

It is important to note that the radiative emission rate near graphene is comparable to 
Fo and therefore negligible compared to SP launching. The carbon sheet is thus eager to 
absorb most of the optical energy released in its vicinity. Consequently, the emitter serves 
as an extremely efficient excitation source of single SPs in graphene. 

A large degree of control over the emission rate can be gained by situating the emitter 
at different distances with respect to the graphene layer, as Fig. 2 illustrates. The near- 
field plots of Fig. 2b, c describe full coupling to SPs at small distance and partial coupling 
at larger separations. The decay and plasmon launching rates both exhibit an exponential 
fall-off with distance predicted by Eq. (Eq. (5)) within the spectral range for which the SPs 
are well defined. At larger energies above the plasmon cutoff, one recovers the same rate as 
in undoped graphene, characterized by a 1/2;^ dependence at large separations (see inset). 

Engineering plasmonic nanostructures 

While we have thus far studied extended graphene sheets, their patterning into nanometer- 
sized cavities yields additional benefits such as extreme field confinement in 3D, engineering 
of resonances, and enhanced coupling efficiency with the far-field. 
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Nanoribbons 

A first example of confinement along one spatial dimension is provided by nanoribbons. The 
results of a thorough theoretical analysis are shown in Fig. 3. These calculations show that 
nanoribbons offer an efficient means of exciting surface SPs in graphene, which can subse- 
quently drive a nearby emitter. Figure 3a depicts the extinction cross section of self-standing 
graphene ribbons for light incident normal to the graphene plane with its polarization across 
the ribbon (supported ribbons lead to similar results, as shown in the Appendix). Plasmon 
confinement is clear from the approximate scaling of the photon energy with the inverse of 
the square root of the ribbon width. The cross section is quite high, demonstrating very 
efficient excitation of SPs. Interestingly, the cross section exceeds the graphene area in some 
cases (black regions). Simultaneously, due to the large wave vector mismatch between SPs 
and far-field, the scattering of plasmons back into photons is very weak, quantified by an 
elastic-scattering contribution to the cross-section that turns out to be more than two orders 
of magnitude smaller than the total cross-section. The latter is therefore dominated by the 
excitation of SPs that are ultimately dissipated in the carbon structure. 

The SP modes exhibit intense field focusing near the edges. The number of induced- 
charge nodes coincides with the order of the mode (1-3), so that odd modes display a net 
dipole moment (e.g., modes 1 and 3 in Fig. 3d). In contrast, even modes (e.g., mode 2) 
couple rather inefficiently to external light because they have higher multipolar character. 

When a line emitter is placed right above the center of the ribbon, odd modes can be 
excited with polarization parallel and across the ribbon (Fig. 3b), whereas even modes couple 
to perpendicular polarization (Fig. 3c). The decay rate is comparable in magnitude to that 
in front of homogeneous graphene (cf. Figures 1 and 3). However, in contrast to homogenous 
graphene, ribbons provide a very efficient way to drive the emitter by external illumination 
because confined plasmons can be excited by an incident plane wave. Effectively the cross 
section at the position of the emitter in Fig. 3b is increased by a factor ^ 750. This reflects 
the ratio of the electric field intensity at the position of the emitter, situated 10 nm away 
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from the carbon sheet, to the incident intensity. 



Nanodisks 

Although ribbons offer an efficient way of exciting SPs, the decay rate of a nearby emitter 
is comparable to a graphene homogeneous sheet, and low photon energies are required. 
In order to boost light-matter interactions at higher photon energies, confinement in all 
dimensions is desirable. Plasmon confinement in all directions is achieved using the circular 
disk cavities illustrated in Fig. 4. This leads to narrow resonances compared to ribbons and 
homogeneous graphene. For simplicity we consider self-standing disks, although like in the 
ribbons, supported disks lead to similar qualitative conclusions. The decay rate (Fig. 4a) is 
significantly boosted at resonance frequencies that can extend up to above Ep. This allows 
one to reach the near-infrared (NIR) region with attainable levels of doping. ^^'^^ 

The light extinction cross-section is also peaked at the SP resonances (Fig. 4b), reaching 
values up to one order of magnitude larger than the disk area. However, similar to the other 
geometries, the radiative emission rate (dashed curves in Fig. 4a) is systematically below 
1% of the total decay rate (sohd curves) in all cases under consideration. This indicates 
that SPs can be efficiently excited by external illumination, but once produced they stay in 
the graphene for up to a few hundred optical cycles (see Q factors below) with negligible 
out-coupling to far-field radiation. 

The modes emerging in the spectra of Fig. 4a have either m = or m = 1 azimuthal 
symmetry when the emitter is situated above the center of the disk and is polarized perpen- 
dicular or parallel to the carbon sheet, respectively. The near electric-field of these modes, 
represented in Fig. 4c, d, clearly shows that the m — 1 plasmon is dipolar, and thus couples 
efficiently to incident light, in contrast to the m — Q plasmon. 

Detailed inspection of the Ep and disk-size dependence of these SPs reveals the following 
properties (see Appendix): the scaling of the plasmon frequency is inherited from the ujp (x 



Ep/Xsp scaling in homogeneous graphene, so that it increases with ^/EF and decreases with 
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the inverse of the square root of the diameter D (in particular, Agp ^ 3D for m = 0, 1); 
maximum Purcell factors F/Fo ^ 10^ — 10^ are consistently obtained; the quality factors 
(extracted from the peak frequency divided by the spectral FWHM) qualitatively follow the 
relation Q ^ uOpT and reach values above 100 for r ^ 10~^^s. 

It should be noted that for such large predicted enhancements of decay rate, the pertur- 
bative treatment of light-matter interactions as described by Eq. (Eq. (10)) can break down, 
giving rise to a new regime of quantum behavior. Physically, the emitter cannot exponen- 
tially decay into the SPs at a rate faster than the plasmonic cavity linewidth = oOp/Q. 
As described below, once Fgp > it should be possible to observe a vacuum Rabi split- 
ting, indicating that an emitted SP can be reversibly and coherently re-absorbed by the 
emitter. ^-^'^^ 

A new regime: Plasmonic vacuum Rabi splitting 

The large Q factors and field concentrations in graphene disks are ideally suited to perform 
quantum optics down to the single-photon level. Here we take a simplified model consisting 
of a two-level quantum emitter (e.g., a quantum dot or a molecule) interacting with a near- 
resonant single mode of a graphene disk. Such a system is characterized by the Jaynes- 
Cummings Hamiltonian^^ H — Ho + i/jnt + ^ext, where 



is the non-interacting part, which includes the plasmon mode of energy hcup and its creation 
and annihilation operators a"^ and a (first term), as well as the unperturbed quantum emitter 
(second term), whose excited state of energy hujQ is created by the operator = |1)(0|, 
connecting its ground |0) and excited |1) levels. We introduce in Hq non-Hermitian damping 
terms to account for inelastic decay channels of both the plasmon mode and the excited 
emitter (e.g., through relaxation and radiative emission) with rates k and Fq, respectively. 




This is consistent with a quantum jump formahsm to describe this open quantum system. 
The SP-emitter interaction is contained in 



and i/ext represents the couphng with an external hght field (see Appendix). 

In this model, the SP-emitter coupling g is determined by comparison with the above 
electromagnetic calculations (Eq. (Eq. (10))), which correspond to the low-coupling limit 
[g <^n)] assuming cjq = cj^, the quantum model yields F = V^ + ^^g^^ / n. In our case, Fq <C 



The properties of the Jaynes-Cummings Hamiltonian are well-studied. When the emit- 
ter and cavity are on resonance and the system is in the strong coupling regime [g > Fq), 
the single-excitation dressed eigenstates consist of even and odd superpositions of the excited 
state of the emitter and a single photon, which can be distinctly resolved. In this regime, an 
initially excited emitter will undergo damped Rabi oscillations at a rate g^ where the emitted 
photon can be re-absorbed before it leaves the cavity. 

The ratio g/i^ reaches a maximum value ^ 4 for 100 nm self-standing graphene disks, 
assuming a reasonable value of the natural decay rate Fq = 5 x 10^ s~^, as we show in Fig. 
5a, b. Actually, g/ ni^ clearly above one for a wide range of doping and disk-size parameters, 
which indicates that the strong-coupling regime is robust. 

A simple signature of strong coupling can be observed in avoided crossings of the ex- 
tinction cross-section of the combined disk-emitter system, a^^^(cj) = Anko lm{a{uj)}^ where 
a{ijj) is the polarizability of the combined system (see Appendix). We show in Fig. 5c the 
extinction cross-section for an emitter of excited energy /icjq = 0.3 eV and natural decay rate 
Fq = 5 X 10^ s"-"^ situated 10 nm above the center of a self-standing 100 nm graphene disk. A 
pronounced vacuum Rabi splitting is observed that can be probed within a single device by 
changing the doping level. 




(6) 
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Conclusion and outlook 

Here, we have described powerful and versatile building blocks for advanced graphene plas- 
monic circuits. These ideas take advantage of the unique combination of extreme field con- 
finement, device tunability and patterning, and low losses that emerge from the remarkable 
structure of graphene and current experimental capabilities for fabrication. These advances 
are expected to both remove a number of obstacles facing traditional metal plasmonics, 
and facilitate new possibilities for manipulating light-matter interactions at the nanoscale 
down to the single-SP level. The simultaneous large bandwidths and field enhancements, 
for example, should enable novel low-power, ultrafast classical or quantum optical devices. 
The strong coupling between single emitters and single SPs could be used to construct 
fast quantum networks, or simulate exotic strongly-interacting condensed matter systems. 
Our proposed techniques could also be potentially applied to manipulate or couple together 
more exotic excitations in graphene, such as thermoplasmons in undoped graphene or 5- 
polarized plasmon modes. Finally, while we have demonstrated the feasibility of graphene 
plasmonics via free-space excitation and detection, the possible applications should be even 
further enhanced with the advent of novel devices such as SP sources, detectors, lasers, opti- 
cal switches and interconnects, plasmon-enhanced photodetectors, and other graphene-based 
nano-optical elements. 
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Figure 1: Coupling of a dipole emitter to doped homogeneous-graphene plasmons. 

a, Near electric field produced by a perpendicular dipole situated 10 nm away from doped 
graphene. The photon and Fermi energies {huj and Ejp) are both 0.5 eV. The real (imagi- 
nary) part of the perpendicular electric field is shown as a red (blue) 3D contour, b, Optical 
dispersion diagram showing the surface plasmon (SP) mode for Ep = 0.5 eV, as well as intra- 
and interband transitions in graphene. c, Real (solid curves) and imaginary (dashed curves) 
parts of the conductivity of doped graphene. d, Decay rate of an excited emitter in front 
of doped graphene as a function of photon emission energy Ep. The rate F is normalized 
to the free-space value Fq. The emission dipole is perpendicular to the graphene and placed 
10 nm away from it. Solid curves show the total decay rate, whereas dashed curves stand 
for the contribution of SP excitation, e, Plasmon dispersion relation in doped graphene. 
The contour plot shows the Fresnel refiection coefficient \rp\ for various values of Ep. The 
dashed lines correspond to the Drude model (Eq. (Eq. (2))). The SP wave vector /cgp ex- 
hibits a quadratic dependence on plasmon energy. The inset shows the propagation distance 
l/Im{/csp} in units of the SP wavelength Agp = 27r/Re{/csp}. The graphene is considered to 
lie on an e = 2 substrate in all cases. 
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Figure 2: Distance dependence of emitter-plasmon coupling, a, Variation of the 
decay rate with distance to the graphene for an emitter polarized perpendicular to the 
carbon sheet. The rate is normalized to the free-space value. Solid (dashed) curves show the 
total (SP-mediated) decay rate. b,c, Near electric-field intensity for two different graphene- 
emitter separations and a Fermi energy Ejp = leV. Foynting vector lines are superimposed 
to the contour plots, with their strength shown in gray scale. The photon energy is 0.5 eV 
and the substrate has e = 2 in all cases. 
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Figure 3: Resonant coupling to graphene ribbons, a, Extinction cross-section of doped 
self-standing graphene ribbons as a function of ribbon width and photon energy for a Fermi 
energy Ejr = 0.2 eV. The hght is incident as shown in the inset. The cross section is normahzed 
to the carbon sheet area. b,c, Decay rate normahzed to free space under the same conditions 
as in a for hue emitters situated 10 nm above the center of the ribbon and polarized parallel 
(b) or perpendicular (c) to it. d, Near electric-field intensity and field lines for the modes 
corresponding to labels 1-3 in a,c. 
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Figure 4: Plasmons in graphene nanodisks. a, Decay rate of an emitter situated 10 nm 
above the center of a doped graphene disk for two different disk diameters D and various 
values of the Fermi energy Ep^ as shown in the legend. Total decay rates (solid curves) are 
compared to the contribution of radiation emission (dashed curves). The emitter is polarized 
parallel to the disk, so that it excites SPs of m = 1 azimuthal symmetry, except in the 
dotted curve, corresponding to perpendicular orientation (m = symmetry) for D = 100 nm 
and Ejp = 0.4 eV. b, Normal-incidence extinction cross-section of the same disks as in a. 
c,d, Near-electric-field intensity of the lowest-energy modes in a = 100 nm disk doped 
to Ep 0.4 eV for m = {hou O.lGeV) and m = 1 {huo 0.26 eV) azimuthal symmetries. 
Electric field-lines are superimposed on the intensity plot. 
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Figure 5: Strong-coupling and vacuum Rabi splitting in graphene nanodisks. a,b, 

Fermi-energy and disk-radius dependence of the strong-coupling parameter g/ nioi an emitter 
placed 10 nm above the center of a doped graphene disk for the first-order (solid curves) and 
second-order (dashed curves) modes with either m = (green curves) or m = 1 (red curves) 
azimuthal symmetries (the two lowest-order modes are shown for each symmetry). The 
natural decay rate of the emitter is Fq = 5 x 10^ s~-^ (see Eq. (Eq. (6)) and below), c, Fermi- 
and photon-energy dependence of the extinction cross-section of a combined emitter-nanodisk 
system under the same conditions as in a. The emitter has a resonance at huo^ = 0.3 eV and 
is oriented parallel to the disk. The cross section is normalized to the maximum resonant 
extinction of the isolated disk. 
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Appendix for 

"Graphene Plasmonics: A Platform for Strong Light- 
Matter Interaction" 

Simulations and relaxation time 

Electromagnetic simulations. We describe ribbons by expanding the induced current 
in Fourier series, assuming that the parallel external field is directed across the ribbon, 
using a supercell with sufl&ciently spaced carbon sheets (see Appendix). Good convergence 
is obtained with 400 Fourier components and supercell spacings of 4 ribbon widths. This 
method produces excellent agreement with an alternative approach fully relying on numerical 
electromagnetic simulations consisting in modeling the graphene as a thin film of dielectric 
function l + Ama/cut and thickness t — 0.5nm, with the edges rounded by hemi-circular 
profiles, for which we find converged results using the boundary element method (BEM), as 
shown in the Appendix. Nanodisks are also simulated with the BEM. The conductivity is 
taken from the fcy ^ hmit of the RPA^^ in all cases. The decay rates are obtained from 
the self-induced field of a dipole using Eq.(3). 

Relaxation time. This is an important parameter because the actual value of r affects 
the plasmon propagation distance, which is nearly proportional to 1/r, although the decay 
rates computed here are rather insensitive to r (because they are A; || -integrated quantities), 
except near the intraband onset (see Appendix). We estimate r from the measured, impurity- 
limited DC mobility i^'i^ /i ^ 10,000cmVVs, which yields r = iiEf/cvI^ ^ lO'^^s for Ef = 
0.1 eV, to be compared with ~ 10~-^^s in gold. We note that this is a very conservative 
value compared to recent observations in high-quality suspended graphene^-^ (/i > 100,000) 
and graphene on boron nitride (/i = 60,000). Optical phonons are known to contribute to 
T above the phonon frequency ^ 0.2 eV. Careful analysis reveals that their effect can be 
incorporated through an effective r ^ 0.5 x 10""^^ s, which produces a reduction in the peak 
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decay rates comparable to the increase in 1/r (e.g., by a factor of 5-10 in the spectra of Fig. 
4, see Appendix). 



Graphene conductivity in the random-phase approximation (RPA) 

The nonlocal conductivity is related to the susceptibility through 

We introduce a finite relaxation time r using the prescription given by Mermin, ^^'^^ which 
preserves the number of charge carriers: 

{l + i/uJT)x{k\\,uj + i/T) 
"'"^^ l + (z/a;rMA;||,a; + z/r)/x(A;||,0)' 



where 



CO 



^ o o/ 



s,s'=± 



. k;r(k|i+k;i) ' 

A;f,|kii+k(,| 



0F{,s'vF\k\\+^\)-9F{svFK) 



u) - 



-VF (sfcJi -s'|k|| +kj||) +iO+ 



is the hnear RPA response function^^'^^ and 9f{E) is the Fermi-Dirac distribution. 

The RPA response admits an analytical expression at zero temperature^^ (i.e., for 6f{E) 
e{EF-E)): 



8kF ^ G(-A-)^[-Re{A_}- 1] + [G(A_) + iw] 6 [Re{ A_} + 1] - G'(A+) 



VfH 



where 



G{z) = log (^z + \Jz'^-l^ 



and A± = {uj / vf ^'2^kF) / k^y Here, the square roots are chosen to yield positive real parts, 
while the imaginary part of the logarithm is taken in the (— 7r,7r] range. Additionally, we 
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have 



e{i-x) 



x\ 1 —x^ 



■COS 



where x = 2kp/k^^. We use these formulas to compute the nonlocal RPA in this document. 

To a good approximation (see below) the conductivity can be evaluated within the local 
RPA (i.e., the fcy limit), which leads to an analytical solution including the dependence 
on T:34 



a{Lj) = 



nh oo + ir 
^2 



+ 



4h 



J log [2cosh{EF/2kBT)] 

4iuj , Hie)-H{ujl2) 

' ^ — 

uj^ — 4e^ 



(7) 



H{uj/2) + —[ de 

TT Jo 



where 



H{s) 



8inh{h£/kBT) 



cosh{EF/kBT) + cosh{hs/kBT) ' 

The first term in Eq. (Eq. (7)) corresponds to intra-band transitions, in which the relaxation 
time has been introduced to make it converge to the Drude model at T = 0. We show below 
that nonlocal effects produce qualitatively similar results as a finite relaxation time r and 
temperature T. Because the actual value of r depends on the quality of the synthetized 
graphene, it can actually be regarded as an effective parameter. Actually, the decay rate of 
an emitter in the vicinity of homogeneous graphene is rather insensitive to the inclusion of 
nonlocal effects and the actual value of r and T within the wide spectral region for which the 
plasmons are well defined (see below). Therefore, we use Eq. (Eq. (7)) for the conductivity in 
all calculations presented here and in the main paper (unless it is stated otherwise) because 
it gives a reasonable description and it is local, so that we assume that it can be also used 
for patterned graphene. 

When T = 0, Eq. (Eq. (7)) reduces to 



a{uj) 



+ 



9{huj-2EF) + -log 

TT 



huj - 2Ef 



huj + 2Ef 



(8) 
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which shows a sudden increase in losses (step function affecting the real part of a) at the 
onset of vertical intra-band transitions, huj = 2Ep. Full inclusion of finite temperature and 
damping leads to a smoother onset, but Eq. (Eq. (8)) contains the main features of the 
graphene conductivity. 




Plasmon energy (eV) 

Figure 6: Spectral dependence of the in-plane plasmon propagation distance (in units of the 
plasmon wavelength) obtained from the RPA for various relaxation times in homogeneous 
graphene supported on an e = 2 material and doped to Ep = leV. 



Fresnel coefficients and plasmon dispersion in homogeneous graphene 

The response of homogeneous graphene is expressed in terms of its Fresnel reflection coef- 
ficients,^^ which can be obtained by applying the customary boundary conditions (AEy = 
AH^ = 0, A(eE^) = 47rcrV • E||/icj, and n x AE^ = (47ra/c)E||) for the fields of incoming p- 
and 5-polarized plane waves as 

ek± — k'\ +4:7rak±k\ luj 



r. = 



ekj_ + fc^ + 4:7rakj_ k'^/uo^ 
k_\_ — k^j_ + Anako / c 
k± + k^j_+ Anako / c ' 
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where e is the permittivity of the substrate on which the graphene is deposited, k{) — uj/c\^ 
the free-space hght wave vector, k±_ = y^/cQ — /cy and k'^ = ^e/cg — A;| are the perpendicular 
wave vectors outside and inside the substrate, respectively, and /cy is the parallel wave vector. 

The dispersion relation of p-polarized surface plasmons (SPs) is determined by the pole 
of Vp^ which yields the equation 

e/^ekl-kl^ + l/^kl-kl^ = -Attcj/uj 

for the plasmon wave vector fcgp. Here, we can use the electrostatic limit of this expression, 

/cgp ~ i{e + l)(jj/ATTa^ 

under the common condition fco <^ l^spl (actually, this condition is fulfilled for huj ^ aEp^ 
where a ^ 1/137 is the fine-structure constant). 

We plot the plasmon dispersion relation in Fig. 1 of the main paper, and we supplement 
it here by showing the 1/e-amplitude-decay propagation length computed from l/Imjfcsp} as 
a function of plasmon energy for various values of r. For the relaxation times considered in 
this work, the susceptibility has an almost linear dependence on the damping rate l/r that 
translates into a linear variation of the plasmon propagation length with this parameter, as 
shown in Fig. Figure 6 for Ep — leV. 

Decay rate and its dependence on conductivity model, tempera- 
ture, and relaxation time 

The decay rate F can be related to the electric field induced by a dipole d on itself E^^^ as^^ 

r = ro + ^im{d*-Ei'^'i}, (10) 
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Figure 7: Decay rate calculated from the RPA at T = for r = 50 fs (solid curves), r = 500 fs 
(crosses), and r = 5,000fs (dots). The emitter is placed 10 nm away from a homogeneous 
graphene sheet deposited on the surface of an e = 2 material. The rate F is normalized to 
the free-space rate Fq. 
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Figure 8: Decay rate obtained from the RPA (symbols) and the local RPA (curves) under 
the same conditions as in Fig. Figure 7 with T = 0, r = jiEp/evp^ and a mobility ji = 
10,000cmVVs (e.g., r ^ lO'^^s for Ep = O.leV). 
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Figure 9: Distance dependence of the decay rate under the conditions of Fig. Figure 8 for a 
photon energy of 0.5 eV. Solid curves: local RPA. Dashed curves: RPA. 
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Figure 10: Decay rate calculated from the local RPA at T = (solid curves), T = 100 
(crosses), and T = 300 (dots) for r = 500 fs under the same conditions as in Fig. Figure 7. 

where Fq = 4/cQ|dp/3/i is the free-space decay rate. When the emitter is above a substrate 
covered with a homogeneous graphene layer, the induced or reflected fleld can be in turn 
related to the Fresnel coeflicients of graphene to yield 



1 

F = Fo + - k\\ dk\\ Re < 



\d\\\\kt^rs-kirp) + 2\d^\^klr. 



^2ik_\_z 



(11) 



where z is the emitter-graphene separation, dy and are the components of the transition 
dipole parallel and perpendicular to the carbon plane, and the integral is extended over 
parallel wave vectors In this work, we use Eq. (Eq. (11)) to compute the decay rate 
in homogeneous graphene, and Eq. (Eq. (10)) for nanoribbons and nanodisks, with E^^^ 
calculated as explained in Sec. ?? and ?? . 

The decay rate F is a -integrated quantity, and therefore, we expect a mild dependence 
on the relaxation time r, except in the neighborhood of the onset of vertical inter-band 
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transitions, where F takes small values that are incremented by the smearing of the electron- 
hole pair (e-h) continuum due to relaxation. This is the conclusion that can be extracted 
from Fig. Figure 7, in which we plot the spectral dependence of F for various values of the 
Fermi energy and we consider a wide range of relaxation parameters. The rate is nearly 
independent of r over the region of existence of surface plasmons and also above the noted 
onset, where it convergences to the undoped graphene level. 

A similar conclusion can be drawn when comparing the RPA with the local RPA (Fig. 
Figure 8). They produce nearly the same results, except in the spectral region extending from 
the plasmon cutoff to the vertical inter-band transition threshold. Clearly, the agreement in 
the plasmon region improves for higher Ejp^ presumably as a result of the momentum cutoff 
for finite separation between the emitter and the graphene, which effectively reduces the 
effect of non-vertical e-h transitions. The two models also differ in the low-energy region. 
The agreement in the plasmonic region is also observed in the distance dependence of F 
(Fig. Figure 9), although severe discrepancies are observed for Ep slightly below the photon 
energy. 

We explore the variation with temperature in Fig. Figure 10. The effect of a finite 
temperature is similar to that of finite relaxation, essentially consisting in smearing the dip 
in the decay rate below the 2Ef onset. 

Fourier expansion method for nanoribbons 

We consider a nanoribbon contained in the z — plane and having translational invariance 
along y. Furthermore, we assume the external field to be independent of y. This is the case 
for illumination with a plane wave normal to the graphene and polarized along x (actually, 
this is the geometry for which we calculate the cross section in this work), and also for 
emission from a line dipole polarized along either x or z and consisting of a continuous 
distribution of identical point dipoles distributed along a line parallel to y (we obtain decay 
rates for this configuration). 
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Under these conditions, the component of the electric field parallel to the graphene is 
along X, and thus, the induced current r]{x) is also along x. The field produced by each 
surface element dxdy is the same as that of a dipole (irj/uj) dxdy. Summing all of these dipole 
contributions, and including the effect of a substrate through its Fresnel coefficients,^^ we 
find the self-consistent relation 

r]{x)/a = El''\x) -^Jdq k^{l-rp) J dxV^(^ -^^//(xO, (12) 

where the x integral is extended over the graphene width, kj_ — ^/cq — g^, ko is the free-space 
light wave vector, and Tp is the Fresnel coefficient of the substrate for p polarization. More 
precisely, Tp — {ek± — k^j_)/{ek± + k^j_)^ where k^j_ — ek'^ — q^. Here, is the external 

electric field along x, which already includes the refiection by the homogeneous dielectric 
substrate of permittivity e. 

We solve Eq. (Eq. (12)) by considering a supercell with the graphene occupying the z = 
and < X < 6 region, and by periodically repeating this unit cell with period a along x. Then, 
we expand the conductivity, the external field, and the surface current in Fourier series. For 
example, the conductivity becomes 

<j{x) = Y.ane'a-^ 

n 

(it is zero outside the graphene and given by Eq. (Eq. (7)) in the graphene), where Qn — 
27rn/a, and the coefficients cr^ are easily obtained from the expansion of the step function 
representing the ribbon. This allows us to project Eq. (Eq. (12)) into Fourier components as 

Vn = - I dxa{x) Ef\x) e-^anx _ _ ^ _ )^^_^,^^,, (13) 

n' 

where /cj and are the same as k±_ and Vp for q = Qn'- Finally, we solve Eq. (Eq. (13)) by 
using standard linear algebra with a finite number of waves M. 
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The scattered near-field is given in terms of the rjn coefficients as 



2tt 



e-*^^-(l + r^)(fc'l,0,5n)A, z<0, 

where the substrate is taken to occupy the z <0 region. 

We obtain convergence for isolated ribbons by taking a equal to 2-4 times b. Then, we 
derive the single-ribbon far-field from the induced current of the ribbon in the first unit cell: 



■gscat 



where ip is the unit vector for scattered p-polarized light, 



27rA;o 



(14) 



is the field amplitude, 6 is the angle relative to the outwards surface normal, — ko {k' — 
ko^/e) above (inside) the substrate, the upper (lower) sign in (Eq. (14)) applies outside 
(inside) the substrate, and 



9n-q 

Here, the reflection coeflicient Vp and the parallel wave vector q — fcosin^ are determined by 
the outgoing angle 9. 

Applying the optical theorem to the transmitted and reflected light upon plane wave 
illumination and using the above expressions for the far fleld, we flnd the extinction cross 
section 



Finally, the decay rate per unit length along y for a line dipole of per-unit-length strength 
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d (_L situated above the substrate is calculated from 

it It 

where the induced field is evaluated at the position of the dipole. 

Convergence of Fourier expansion and boundary element method 
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Figure 11: Convergence of the Fourier expansion (FE) and the boundary element method 
(BEM) for a graphene nanoribbon of width 100 nm and Ep = 0.2 eV. We represent the extinc- 
tion cross section calculated with various values of the convergence parameters for normally- 
incident light. The number of Fourier components M and boundary parametrization points 
are given in the text labels. 

The Fourier expansion method converges when the number of Fourier coefficients is in- 
creased, as we show in Fig. Figure 11 (symbols). This method produces results in excellent 
agreement with an alternative approach fully relying on numerical simulations, consisting 
in modeling the graphene as a thin film of dielectric function l + Ania/ujt and thickness t, 
with the edges rounded by hemi-circular profiles, for which we find converged electromag- 
netic results using the boundary element method (BEM),^^ as shown in Fig. Figure 11 (solid 
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curves). 

The agreement between the semi-analytical Fourier expansion and the BEM confirms the 
validity of the latter to describe graphene as a thin effective layer of dielectric. Actually, 
we use this method to simulate graphene disks, because an analytical expansion for them 
becomes too involved and does not add much insight into the problem. 

Supported vs self-standing nanoribbons 
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Figure 12: (a) Extinction cross section of doped graphene ribbons deposited on an e = 2 
material as a function of ribbon width and photon energy for a Fermi energy Ejr = 0.2 eV. 
The hght is incident as shown in the inset. The cross section is normahzed to the carbon 
sheet area. (b,c) Decay rate normahzed to free space under the same conditions as in (a) 
for a hne emitter situated 10 nm above the center of the ribbon and polarized either parallel 
(b) or perpendicular (c) to it. 



We show in Fig. Figure 12 calculations similar to those of Fig. 3 of the main paper, but 
for graphene ribbons supported on the surface of an e = 2 material rather than self-standing. 
The results are qualitatively the same for supported and for self-standing graphene. The 
extinction cross section and the decay rates have similar magnitude in both cases. The 
only difference that is worth noticing is the redshift in the plasmon energy in the supported 
graphene, which is consistent with the scaling of oc 1/ ^/e + l predicted by the Drude model. 
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Distance dependence of the decay rate near a nanodisk 

In the main paper, we have discussed the decay rate for an emitter situated at a fixed point 
along the axis of self-standing circular graphene disks. Here we consider the variation of the 
peak decay rate as a function of position of the emitter. Figure Figure 13a shows the variation 
of the rate with distance to the graphene along the axis of the disk. The rate decays with 
distance z roughly as exp(— 47rz/Asp) (i.e., it follows the same exponential attenuation as in 
homogeneous graphene), where the plasmon wavelength Agp is 290 nm for the m = 1 mode 
and 94 nm for the m = mode. The variation along parallel displacements (Fig. Figure 13b) 
is less trivial, but it qualitatively follows the near-field intensities shown in Fig. 4c, d of the 
main paper. 
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Figure 13: Position dependence of the decay rate F near a self-standing graphene disk of 
diameter 100 nm and Ejp = 0.4 eV. The rate is given as a function of position along the dashed- 
line excursions shown in the upper insets. Plot (b) is obtained for a distance of 10 nm from 
the graphene plane. The rate is calculated at the peak position of the lowest-order m = 1 
and m = resonances, respectively (see Fig. 4 of the main paper). 
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Figure 14: Plasmons in graphene nanodisks: Size and doping dependence. (a,b) 

Fermi-energy and disk-radius dependence of the SP energy in doped graphene disks for 
the first- (sohd curves) and second-order (dashed curves) modes with either m = (green 
curves) or m = 1 (red curves) azimuthal symmetries. The two lowest-order modes are shown 
for each symmetry. (c,d) Decay rate F of an emitter located 10 nm above the center of the 
disk, normalized to the rate in free space Fq, and calculated at the energies of the SPs in 
(a,b). The emitter is polarized parallel (perpendicular) to the disk in the red (green) curves. 
(e,f) Quality factor of the resonances considered in (a,b). 
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Figure 15: Plasmons in graphene nanodisks: Relaxation dependence. Decay rate of 
an emitter situated 10 nm above the center of a doped graphene disk for two different disk 
diameters D and two values of the Fermi energy as shown by text labels. The emitter 
is polarized parallel to the disk. Dashed curves are calculated with r ~ 0.5 X 10-1^ s. Solid 
curves are obtained with r = 4 x 10~^^s for Ep = 0.4 eV and r = 6 x 10~^^s for Ep = 0.6 eV. 

Plasmons in graphene nanodisks: Size, doping, and relaxation de- 
pendence 

Figure Figure 14 summarizes the Ep and disk-size dependence of SPs in graphene nanodisks. 
The evolution of the plasmon energy is inherited from the ujp oc Ep/\^^ scaling in homo- 
geneous graphene, so that it increases with ^/EF (Fig. Figure 14a) and decreases with the 
inverse of the square root of the radius (Fig. Figure 14b). 

The peak decay rate has a weak dependence on both Ep and disk size (Fig. Figure 14c, d), 
essentially reflecting the increase in SP lifetime when its energy is positioned close to Ep. 
Maximum rates ~ 10^ Fq are consistently obtained near these disks at the distance z = 10 nm 
chosen in the calculations. 

The quality factor of the SP resonances Q (Fig. Figure 14e,f, extracted from the peak 
frequency divided by the FWHM) shows a strong increase with Ep^ in agreement with 
the longer propagation distance observed in homogeneous graphene at higher doping. Our 
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reported values of Q > 100 are consistent with the moderate relaxation time employed in the 
calculations, r ^ 10"-^^ s (i.e., Q ^ oopr). Actually, the values of Q reported in Fig. Figure 14e,f 
are qualitatively well predicted by this formula using the energies of Fig. Figure 14a,b as 
input. At the same time, much higher values of the mobility fi have been reported, which 
should lead to larger Q's and peak rates in direct proportion to /i, at least below the optical 
phonon frequency. 

The role of plasmons in the relaxation of graphene is not yet well understood, although 
careful analysis reveals that their effect can be incorporated through an effective r ^ 0.5 x 
10~-^^s. We show in Fig. Figure 15 results for the decay rate near nanodisks obtained with 
this effective value of r (dashed curves), compared to calculations obtained with r = jiEp/evp 
and a mobility ji = 10,000cm^/Vs (solid curves, taken from Fig. 4 of the main paper). The 
shorter relaxation time due to phonons produces a decrease in both the peak rate and Q by 
a factor of ^5-9 in the spectra shown here, which increases linearly with Ep and decreases 
smoothly with uj. 

Field lines in near-field plots 

The near-electric-fields shown in the main paper for ribbons and disks are obtained close to 
a non-degenerate resonance at = cjq, so that they take the form E^^^ + F/(cJo ~ ^ ~ '^7/2), 
where 7 is the plasmon relaxation rate. Since the ribbon width is much smaller than the 
light wavelength, the near field must be almost electrostatic, and therefore, E and F are 
approximately real vectors. Interestingly, the on-resonance induced field z2F/7 is almost 
imaginary, in good agreement with our numerical simulations. Using a real transition dipole, 
the induced fields plotted in the main paper are more than 99% imaginary, and this actually 
allows us to extract the field fines that are shown there (i.e., field fines corresponding to a 
nearly real, electrostatic electric field). 
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Polarizability of the combined SP-emitter system in the Jaynes- 
Cummings model 

The interaction with the external field E{t) can be written 

which involves the plasmon and emitter dipole operators Pp — dpa + d^a'^ and Pq = doa + 
dQa^^ where dp and do are their respective transition dipoles. 

We derive the polarizability a from the first-order-perturbation-theory steady-state so- 
lution of the model Hamiltonian for a faint external field of the form E{t) — 2Re{£'oe~^^^}. 
The polarizability is defined through the relation p{t) = 2Re{a{uj)Eoe~^^^}^ obtained from 
the expected value of the combined induced dipole, p = {Pp + Pq). After some algebra, we 
find 

where 

_ {uQ-u- iTo/2)\dp\^ + {ujp-uj- m/2) |q!q ^ - 2glm{dpd'^} 
(cjQ — uj — iVo / 2) {ujp — uj — iti/2) — g'^ 

In the calculations presented in Fig. 5c of the main paper we assume do to be negligible 
compared to dp. Furthermore, the model parameters cj^, F, and Q are extracted from 
the second m = 1 mode of Fig. Figure 14a, c,e (red dashed curves) following the procedure 
described in the main paper. 
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